This script was most recently run on Fri Feb 12 14:53:12 2016.
combined_data <- rbind(cohort7_filter, cohort8_filter)
combined_data$Day <- as.factor(combined_data$Day)
combined_data$Interval <- as.factor(combined_data$Interval)
combined_data$Light_Dark <- as.factor(combined_data$Light_Dark)
ggplot(data=combined_data[combined_data$Animal_ID=="1767",], aes(x=Date.Time, y=Weight, group=Animal_ID, colour=Animal_ID)) +
geom_point() + geom_line()+
xlab("Interval") +
ylab("Food (g)")
ggplot(data=combined_data[combined_data$Animal_ID=="2112",], aes(x=Date.Time, y=Weight, group=Animal_ID, colour=Animal_ID)) +
geom_point() + geom_line()+
xlab("Interval") +
ylab("Food (g)")
ggplot(data=combined_data[combined_data$Animal_ID=="2105",], aes(x=Date.Time, y=Weight, group=Animal_ID, colour=Animal_ID)) +
geom_point() + geom_line()+
xlab("Interval") +
ylab("Food (g)")
ggplot(data=combined_data[combined_data$Animal_ID=="2106",], aes(x=Date.Time, y=Weight, group=Animal_ID, colour=Animal_ID)) +
geom_point() + geom_line()+
xlab("Interval") +
ylab("Food (g)")
ggplot(data=combined_data[combined_data$Animal_ID=="2118",], aes(x=Date.Time, y=Weight, group=Animal_ID, colour=Animal_ID)) +
geom_point() + geom_line()+
xlab("Interval") +
ylab("Food (g)")
ggplot(data=combined_data[combined_data$Animal_ID=="1757",], aes(x=Date.Time, y=Weight, group=Animal_ID, colour=Animal_ID)) +
geom_point() + geom_line()+
xlab("Interval") +
ylab("Food (g)")
ggplot(data=combined_data[combined_data$Animal_ID=="2119",], aes(x=Date.Time, y=Weight, group=Animal_ID, colour=Animal_ID)) +
geom_point() + geom_line()+
xlab("Interval") +
ylab("Food (g)")
ggplot(data=combined_data[combined_data$Animal_ID=="1753",], aes(x=Date.Time, y=Weight, group=Animal_ID, colour=Animal_ID)) +
geom_point() + geom_line()+
xlab("Interval") +
ylab("Food (g)")
ggplot(data=combined_data[combined_data$Animal_ID=="1760",], aes(x=Date.Time, y=Weight, group=Animal_ID, colour=Animal_ID)) +
geom_point() + geom_line()+
xlab("Interval") +
ylab("Food (g)")
ggplot(data=combined_data[combined_data$Animal_ID=="1754",], aes(x=Date.Time, y=Weight, group=Animal_ID, colour=Animal_ID)) +
geom_point() + geom_line()+
xlab("Interval") +
ylab("Food (g)")
#remove data from animals 2105 and 2119
combined_data <- combined_data[! combined_data$Animal_ID %in% c(2105, 2119, 1753, 1754, 1757, 1760),]
se <- function(x) sd(x)/sqrt(length(x))
combined_summary <-
combined_data %>%
group_by(Cohort, Animal_ID, Light_Dark) %>%
summarise(
sum_weight = sum(Weight),
mean_weight = mean(Weight),
se_weight = se(Weight),
sum_duration = sum(Duration),
mean_duration = mean(Duration),
se_duration = se(Duration),
num_feed = length(Weight),
food_per_sec = sum_weight/sum_duration
#weight.shapiro = shapiro.test(Weight)$p.value
)
combined_summary <- combined_summary[!is.na(combined_summary$Light_Dark), ]
#Number of feedings for each mouse
ggplot(data=combined_summary, aes(x=Animal_ID, y=num_feed, colour=Animal_ID)) +
geom_point(data=combined_summary, aes(fill=Animal_ID)) +
facet_grid(.~Light_Dark) +
ggtitle("Feeding habit")+
xlab("Animal_ID") +
ylab("Number of feedings")
ggsave(file="figure/Cohort7_8_Number_of_feeding.pdf")
#Avg Duration per feeding
ggplot(data=combined_summary, aes(x=Animal_ID, y=mean_duration, colour=Animal_ID)) +
geom_bar(position=position_dodge(), aes(fill=Animal_ID),stat="identity") +
geom_errorbar(aes(ymin=mean_duration-se_duration, ymax=mean_duration+se_duration), width=0.25, position=position_dodge(.9)) +
facet_grid(.~Light_Dark) +
ggtitle("Feeding Duration")+
xlab("Animal_ID") +
ylab("Avg duration per feeding (s)")
ggsave(file="figure/Cohort7_8_duration_per_feeding.pdf")
ggplot(combined_summary, aes(x=Animal_ID, y=sum_weight, group=Animal_ID, colour=Animal_ID)) +
geom_point()+
facet_grid(Light_Dark~.)+
ggtitle("Cohort7 total food intake")+
xlab("Animal ID")+
ylab("Food (g)")
ggsave("figure/Cohort7_8_total_food_intake_in_each_animal.pdf")
#calculate the amount of food intake for each day
combined_each_animal_by_day <-
combined_data %>%
group_by(Cohort, Animal_ID, Day) %>%
summarise(
sum_weight = sum(Weight), #get the food for each day
sum_duration = sum(Duration),
mean_duration = mean(Duration),
se_duration = se(Duration),
num_feed = length(Weight)
)
#Is the amount of food intake per day the same in all animals?
#NO, each animal consume different amount of food per day
kruskal.test(sum_weight ~ Animal_ID, data=combined_each_animal_by_day)
##
## Kruskal-Wallis rank sum test
##
## data: sum_weight by Animal_ID
## Kruskal-Wallis chi-squared = 113.4231, df = 44, p-value =
## 4.826e-08
kruskal.test(sum_duration ~ Animal_ID, data=combined_each_animal_by_day)
##
## Kruskal-Wallis rank sum test
##
## data: sum_duration by Animal_ID
## Kruskal-Wallis chi-squared = 111.3729, df = 44, p-value =
## 9.263e-08
#calculate the amount of food intake for each animal from the each day statistics
combined_each_animal <-
combined_each_animal_by_day %>%
group_by(Cohort, Animal_ID) %>%
summarise(
cum_food = sum(sum_weight),
mean_food = mean(sum_weight),
se_food = se(sum_weight),
num_day = length(Day),
avg_duration = mean(mean_duration),
se_duration = mean(se_duration),
num_feed = sum(num_feed),
food.day = cum_food/num_day
#weight.shapiro = shapiro.test(sum_weight)$p.value
)
#remove animals with cummulative food intake < 1 g
sub_data <- combined_each_animal[combined_each_animal$cum_food>1,]
ggplot(data=sub_data, aes(x=Animal_ID, y=mean_food, colour=Animal_ID)) +
geom_bar(position=position_dodge(), aes(fill=Animal_ID),stat="identity") +
geom_errorbar(aes(ymin=mean_food-se_food, ymax=mean_food+se_food), width=0.25, position=position_dodge(.9)) +
ggtitle("Food intake per day")+
xlab("Animal_ID") +
ylab("Food (g/day)")
ggsave("figure/Food_intake_per_day.pdf")
## Saving 7 x 5 in image
#Number of feedings for each mouse
ggplot(data=sub_data, aes(x=Animal_ID, y=num_feed, colour=Animal_ID)) +
geom_point(data=sub_data, aes(fill=Animal_ID)) +
ggtitle("Feeding habit")+
xlab("Animal_ID") +
ylab("Number of feedings")
ggsave("figure/Number_of_feedings.pdf")
## Saving 7 x 5 in image
#Avg Duration per feeding
ggplot(data=combined_each_animal, aes(x=Animal_ID, y=avg_duration, colour=Animal_ID)) +
geom_bar(position=position_dodge(), aes(fill=Animal_ID),stat="identity") +
geom_errorbar(aes(ymin=avg_duration-se_duration, ymax=avg_duration+se_duration), width=0.25, position=position_dodge(.9)) +
#facet_grid(.~Light_Dark) +
ggtitle("Feeding habit")+
xlab("Animal_ID") +
ylab("Avg duration per feeding per day (s/day)")
ggsave("figure/Cohor7_duration_per_feeding_per_day.pdf")
## Saving 7 x 5 in image
ggplot(sub_data, aes(x=Animal_ID, y=cum_food, group=Animal_ID, colour=Animal_ID)) +
geom_point()+
ggtitle("Total food intake")+
xlab("Animal ID")+
ylab("Food (g)")
ggsave("figure/Total_food_intake_in_each_animal.pdf")
## Saving 7 x 5 in image
food_dark <- combined_summary[combined_summary$Light_Dark=="Dark",]
food_light <- combined_summary[combined_summary$Light_Dark=="Light",]
dark.shapiro <- shapiro.test(food_dark$sum_weight)
wilcox.test(food_dark$sum_weight, food_light$sum_weight, alternative="greater", conf.int=T)
##
## Wilcoxon rank sum test
##
## data: food_dark$sum_weight and food_light$sum_weight
## W = 1222, p-value = 0.04589
## alternative hypothesis: true location shift is greater than 0
## 95 percent confidence interval:
## 0.15 Inf
## sample estimates:
## difference in location
## 4.08
boxplot(sum_weight~Light_Dark, data=combined_summary, ylab="Food (g)")
boxplot(sum_duration~Light_Dark, data=combined_summary, ylab="Duration (s)")
wilcox.test(sum_weight ~ Light_Dark, data=combined_summary, alternative="greater")
##
## Wilcoxon rank sum test
##
## data: sum_weight by Light_Dark
## W = 1222, p-value = 0.04589
## alternative hypothesis: true location shift is greater than 0
wilcox.test(sum_duration~Light_Dark, data=combined_summary, alternative="greater")
##
## Wilcoxon rank sum test
##
## data: sum_duration by Light_Dark
## W = 1266, p-value = 0.02039
## alternative hypothesis: true location shift is greater than 0
library(lattice)
## Warning: package 'lattice' was built under R version 3.1.3
xyplot(sum_weight~Day|Animal_ID, data=combined_each_animal_by_day, xlab = "Day", ylab = "Food (g)")
The results suggested that mice may consumed a bit more food in the dark and light with a border Wilcoxon’s p-value of 0.0458865. Mice spent more time eating in the Dark than in the Light, Wilcoxon’s p-value = 0.0203912.
#take the average food per day of each animal = the average food intake for each cohort
data_per_cohort <-
sub_data %>%
group_by(Cohort) %>%
summarise(
num_animals = length(Animal_ID),
avg_food_per_day = mean(mean_food),
se_food_per_day= sqrt(sum(se_food^2)/num_animals),
avg_duration_cohort = mean(avg_duration),
se_duration_cohort = sqrt(sum(se_duration^2)/num_animals),
num_feed = sum(num_feed)
)
shapiro.test(sub_data[sub_data$Cohort=='7',]$mean_food)
##
## Shapiro-Wilk normality test
##
## data: sub_data[sub_data$Cohort == "7", ]$mean_food
## W = 0.7066, p-value = 4.697e-05
wilcox.test(mean_food ~ Cohort, data=sub_data, alternative="two.sided")
##
## Wilcoxon rank sum test
##
## data: mean_food by Cohort
## W = 170, p-value = 0.06909
## alternative hypothesis: true location shift is not equal to 0
ggplot(data=data_per_cohort, aes(x=Cohort, y=avg_food_per_day, colour=Cohort)) +
geom_bar(position=position_dodge(), aes(fill=Cohort),stat="identity") +
geom_errorbar(aes(ymin=avg_food_per_day-se_food_per_day, ymax=avg_food_per_day+se_food_per_day), width=0.25, position=position_dodge(.9)) +
ggtitle("Food Intake")+
xlab("Cohort") +
ylab("Daily food intake (g)")
ggsave("figure/Daily_food_intake_in_each_cohort.pdf")
## Saving 7 x 5 in image
shapiro.test(sub_data[sub_data$Cohort=='8',]$avg_duration)
##
## Shapiro-Wilk normality test
##
## data: sub_data[sub_data$Cohort == "8", ]$avg_duration
## W = 0.7624, p-value = 5.713e-05
wilcox.test(avg_duration~Cohort, data=sub_data, alternative="two.sided")
##
## Wilcoxon rank sum test
##
## data: avg_duration by Cohort
## W = 280, p-value = 0.5042
## alternative hypothesis: true location shift is not equal to 0
ggplot(data=data_per_cohort, aes(x=Cohort, y=avg_duration_cohort, colour=Cohort)) +
geom_bar(position=position_dodge(), aes(fill=Cohort),stat="identity") +
geom_errorbar(aes(ymin=avg_duration_cohort-se_duration_cohort, ymax=avg_duration_cohort+se_duration_cohort), width=0.25, position=position_dodge(.9)) +
ggtitle("Duration")+
xlab("Cohort") +
ylab("Duration per feeding (s/day)")
ggsave("figure/Avg_duration_per_feeding_per_day_in_each_cohort.pdf")
## Saving 7 x 5 in image
ggplot(data=sub_data, aes(x=Cohort, y=num_feed, colour=Cohort)) +
geom_boxplot() +
geom_jitter() +
xlab("Cohort") +
ylab("Number of feedings")
ggplot(data=sub_data, aes(x=Cohort, y=avg_duration, colour=Cohort)) +
geom_boxplot() +
geom_jitter() +
xlab("Cohort") +
ylab("Duration per feeding per day (s)")
ggplot(data=sub_data, aes(x=Cohort, y=mean_food, colour=Cohort)) +
geom_boxplot() +
geom_jitter() +
xlab("Cohort") +
ylab("Daily food intake (g)")
The food intake in Cohort 7 is not normally distributed as Shapiro-Wilk’s p-value is 4.697276810^{-5}. The duration per feeding per day in cohort 8 is not normally distributed as tested by Shapiro-Wilk test, 5.712689810^{-5}. The amount of daily food intake in cohort 7 and 8 are not different, Wilcoxon p-value is 0.0690876. The average duration per feeding per day are not different between cohort 7 and 8 as tested by Wilcoxon rank sum test (p = 0.5041839).
body_comp_7_file <- "../data/processed/body_composition_cohort7.csv"
body_comp_7 <- read.csv(body_comp_7_file)
body_comp_8_file <- "../data/processed/body_composition_cohort8.csv"
body_comp_8 <- read.csv(body_comp_8_file)
sub_body_comp_8 <- body_comp_8[, !colnames(body_comp_8)=='Grip.Strength..4.Paw.']
combined_body_comp <- rbind(body_comp_7, sub_body_comp_8)
food_LBM <- merge(combined_body_comp, sub_data, by.x="animal.MouseID", by.y="Animal_ID")
sub_food_LBM <- food_LBM[food_LBM$Fasting_stat %in% c('Pre_fasting', 'After_HFD'),]
food_per_LBM <-
sub_food_LBM %>%
group_by(animal.MouseID, Fasting_stat) %>%
summarise(
food.LBM = cum_food/Lean,
food.BW = cum_food/Weight,
food.FM = cum_food/Fat,
food.LBM.day = (mean_food/Lean),
food.BW.day = (mean_food/Weight),
food.FM.day = mean_food/Fat,
food.day = food.day
)
food_per_LBM$animal.MouseID <- as.factor(food_per_LBM$animal.MouseID)
food_per_LBM$Fasting_stat <- relevel(as.factor(food_per_LBM$Fasting_stat), ref='Pre_fasting')
ggplot(data=food_per_LBM, aes(x=animal.MouseID, y=food.BW.day, colour=animal.MouseID)) +
geom_bar(position=position_dodge(), aes(fill=animal.MouseID),stat="identity") +
facet_grid(.~Fasting_stat)+
xlab("Animal_ID") +
ylab("Food intake (g/g body weight per day)")
with(food_LBM, plot(Weight,food.day, pch=19, las=1,
col=c('green','red'),
xlab="Fat Mass (g)",
ylim=c(0,max(food.day)),
ylab="Daily Food Intake (g)"))
cor(food_LBM$Lean, food_LBM$food.day, method="kendall")
## [1] -0.01660486
cor.test(food_LBM$Lean, food_LBM$food.day, method="kendall")
##
## Kendall's rank correlation tau
##
## data: food_LBM$Lean and food_LBM$food.day
## z = -0.2835, p-value = 0.7768
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
## tau
## -0.01660486
HFD <- food_LBM[food_LBM$Fasting_stat=="After_HFD",]
with(HFD, plot(Weight,food.day, pch=19, las=1,
#col=c('green','red','blue'),
xlab="Body Weight (g)",
ylim=c(0,max(food.day)),
ylab="Daily Food Intake (g)"))
cor(HFD$Weight, HFD$food.day, method="kendall")
## [1] 0.03945376
cor.test(HFD$Weight, HFD$food.day, method="kendall")
## Warning in cor.test.default(HFD$Weight, HFD$food.day, method = "kendall"):
## Cannot compute exact p-value with ties
##
## Kendall's rank correlation tau
##
## data: HFD$Weight and HFD$food.day
## z = 0.3816, p-value = 0.7028
## alternative hypothesis: true tau is not equal to 0
## sample estimates:
## tau
## 0.03945376
There is no correlation between the food intake before fasting and the weight after HFD. The correlation is 0.0394538 with the p-value of 0.901549.
## R version 3.1.2 (2014-10-31)
## Platform: x86_64-apple-darwin13.4.0 (64-bit)
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] lattice_0.20-33 lubridate_1.5.0 dplyr_0.4.3 plyr_1.8.3
## [5] ggplot2_2.0.0
##
## loaded via a namespace (and not attached):
## [1] assertthat_0.1 colorspace_1.2-6 DBI_0.3.1 digest_0.6.9
## [5] evaluate_0.8 formatR_1.2.1 grid_3.1.2 gtable_0.1.2
## [9] htmltools_0.3 knitr_1.12.3 labeling_0.3 lazyeval_0.1.10
## [13] magrittr_1.5 munsell_0.4.2 parallel_3.1.2 R6_2.1.2
## [17] Rcpp_0.12.3 reshape2_1.4.1 rmarkdown_0.9.2 scales_0.3.0
## [21] stringi_1.0-1 stringr_1.0.0 tools_3.1.2 yaml_2.1.13